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ABSTRACT 

The oscillation modes of a simple polytropic stellar model are studied. Using a new 
numerical approach (based on integration for complex coordinates) to the problem 
for the stellar exterior we have computed the eigenfrequencies of the highly damped 
w-modes. The results obtained agree well with recent ones of Leins, NoUert and Sof- 
fel [Phys. Rev. D 48 3467 (1993)]. Specifically, we are able to explain why several 
modes in this regime of the complex frequency plane could not be identified within 
the WKB approach of Kokkotas and Schutz [Mon. Not. R. astr. Soc. 255 119 (1992)]. 
Furthermore, we have established that the "kink" that was a prominent feature of 
the spectra of Kokkotas and Schutz, but did not appear in the results of Leins et ai, 
was a numerical artefact. Using our new numerical code we are also able to compute, 
for the first time, several of the slowly damped (p) modes for the considered stellar 
models. For very compact stars we find, somewhat surprisingly, that the damping 
of these modes does not decrease monotonically as one proceeds to higher oscillation 
frequencies. The existence of low-order modes that damp away much faster than antic- 
ipated may have implications for questions regarding stellar stability and the lifetime 
of gravitational-wave sources. The present results illustrate the accuracy and reliabil- 
ity of the complex-coordinate method and indicate that the method could prove to 
be of great use also in problems involving rotating stars. There is no apparent reason 
why the complex-coordinate approach should not extend to rotating stars, whereas it 
is accepted that all previous methods will fail to do so. 
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1 INTRODUCTION 



In the concluding remarks of their pioneering paper on non- 
radial oscillations of neutron stars Thorne and CampoUataro 
(1967) described it as "just a modest introduction to a story 
which promises to be long, complicated and fascinatinrj'' . The 
story has undoubtedly proved to be intriguing, and many 



authors have contributed to our present understanding of 
the pulsations of neutron stars. Of special interest are the 
attempts to calculate the eigenfrequencies of the so-called 
quasinormal modes of a stellar system. These are solutions 
to the perturbation equations which agree with the neces- 
sary boundary conditions at the centre and at the surface 
of the star, and at the same time behave as purely outgoing 
waves at spatial infinity. In general relativity pulsations of 
the stellar fluid will normally be damped due to the emis- 
sion of gravitational radiation. Hence, the characteristic fre- 
quencies of the system will have complex values. The real 
parts correspond to physical oscillation frequencies, while 



the imaginary parts describe the damping rate of each os- 
cillation mode. 

The first attempts to compute outgoing-wave modes for 
relativistic stars concerned modes that have an analogue 
in Newtonian theory. It was believed (and verified by the 
early calculations for nonrotating stars) that the Newtonian 
modes of oscillation are shifted only slightly because of the 
coupling to gravitational radiation: Each characteristic fre- 
quency adopts a very small imaginary part. In practice, it 
is not straightforward to determine these modes with great 
accuracy. This diffic ulty became o bvious with the first nu- 
merical calculations (Thorne 1969). An approach based on a 



variational principle devised by Detweiler and I pser (1973) 
also suffers from this problem (Detweiler 1975). The diffi- 



culty of determining small imaginary parts with some ac- 
curacy was not over come until much later when Lindblom 
and Detweiler ( 1985 ) combined numerical integration of the 
equations for the interior with a numerically integrated so- 
lution for the exterior. Using this technique they computed 
the fundamental (f) mode frequency for many realistic neu- 
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tron star equations of state. However, they later realized that 
the fourth-order system of equations they had used for the 
stellar interior could become singular. Consequently, they 
sugges ted that a different choice o f dependent variables be 
made (Detweiler & Lindblom 1981:). Also w orth m entioning 



in this context is a paper by Balbinski et al. (1985) where re- 



sults obtained by the Lindblom-Detweiler scheme were found 
to be in satisfactory agreement with the predictions of the 
standard quadrupole formula. 

The last few years have seen several important contribu- 
tions to this field. In an impressive series of papers, Chan- 
drasekhar and Ferrari have reformulated the problem. In 
their final syste m of equations all fluid perturba tions have 
been eliminated ( Chandrasekhar & Ferrari 1991a). This has 



the advantage that one can view the problem of gravitational 
waves being scattere d by a star in a way similar to that 
used for black holes (Chandrasekhar 1983). In the second 
paper of the series, Chandrasekhar and Ferrari show that 
odd-parity (axial) modes can be excited by the even-parity 
(polar) modes in the slow rotation limit (Chandrasekhar & 
Ferrari 1991b). The coupling is due to the Lense-Thirring 



effect. In the following papers, they demonstrate that their 
algorithm for computing eigenfrequencies can be trusted so 
long as th e imaginary part is consider ably smaller than the 
real part (Chandrasekhar et al. 1991), that the odd-parity 



(axial) modes can be of importance fo r extremely compact 
stars ( Chandrasekhar fc Ferrari 1991c ) and that the idea of 
a complex ang ular momentum can be useful in astrophysi- 
cal situations ( Chandrasekhar fc Ferrari 19921). A summary 
of their work has been written by Ferrari (1992). Another 
recent reformulation of the equations governing a per turbed 
stellar model was proposed by Ipser and Price (1991). They 
discuss the relation between the Regge- Wheeler gauge (and 
the resulting equations) and the diagonal gauge used by 
Chandrasekhar and Ferrari (Price fc Ipser 1991). In fact, 
they reformulate the stellar perturbation equations in such a 
way that all fluid perturbations can be eliminated also in the 
Regge- Wheeler gauge. An interesting question is whether 
these new formulations of the pulsation problem will make 
computation of mode-frequencies any easier. We will not ap- 



proach that issue here, but feel that an investigation of that 
kind is of great importance and should be encouraged. 

Another chapter of the story begins with t he stu dy of 
a very simple toy model. Kokkotas and Schutz ( 1986 ) sug- 
gested that a plausible model for a stellar system would be a 
finite string - representing the fluid of the star - coupled by 
means of a spring to a semi-infinite string - the dynamical 
spacetime. They found that such systems could accommo- 
date a new family of oscillation modes. These would be asso- 
ciated with the gravitational degrees of freedom, rather than 
the pulsations of the fluid, and would be strongly damped. 
Recently, a sli ghtly more realistic toy model led Baumgarte 
and Schmidt (1993) to much the same c onclus ions. Yet an- 
other toy model, suggested by Kokkotas (1985), shows a lot 
of similarities with the axial (o dd pa rity) modes discussed 



by Chandrasekhar and Ferrari (1991) and has recently led 
t o the discovery of a branch of strongly damped axial modes 
([Kokkotas 1994|)- 



A computation of highly damped modes for a realis- 
tic stellar model is not trivial, however. The main difficulty 
involves numerically separating the ingoing and outgoing- 
wave solutions at spatial infinity: The ingoing solution dies 



exponentially as r ^ cxj while the outgoing one grows. This 
problem is well-known from studies of n uasinormal modes 
for black holes (see Andersson et al. ( |l993[ ) for a detailed dis- 
cussion). It seems reasonable to assume that methods that 
have proved useful in studies of the black-hole problem can 
be adapted to the stellar situation. A series expansion ap- 
proach based on a four-term recurrence r elation - sim ilar to 
that developed by Leaver for black holes ( Leaver 1985[ ) - was 
used by Kojima (1988) to verify that the strongly damped 



modes (referred to as w-modes because of their connection 
to the gravitational waves) indeed do exi st for realistic stel- 
lar models. Kokkotas and Schutz (1992) used a WKB ap- 
proach (essentially a geometrical optics assumption of no 
refiection of waves in the exterior vacuum) in their calcula- 
tions of w-mode spectra for several models. Their main re- 
sults have recently been verified by Leins, NoUert and Soffel 
( [l993| ). ' n their calculations, Leins et al. employed two dif- 



ferent approac hes for the e xterior: Leaver's continued frac- 
tion approach (Leaver 1985) and a Wronskian technique that 
has proved extremely pow erful for Schwarzschild black holes 
( [Nollert fc Schmidt 199^ ). 

Although the results obtained by Leins et al. generally 
agree with those of Kokkotas and Schutz there are some 
differences. Leins et al. found new modes with considerably 
smaller oscillation frequencies and also higher damping than 
those that had been found by the WKB-technique. Their 
calculations also suggest that a "kink" that was apparent 
in the spectra of Kokkotas and Schutz does not exist. At 
first sight, one would be tempted to believe that the cal- 
culations of Leins et al. are the most reliable since the two 
approaches they used to deal with the exterior do not, in 
principle, involve any approximations. This, however, does 
not explain why the WKB approach becomes less reliable as 
the oscillation frequency increases. In fact, this is contrary 
to all expectations: The geometrical optics argument should 
be valid for high frequencies. 

In this paper we attempt to settle these matters. We 
will combine a numerical phas e- amp litude approach (simi- 
lar to that used by Andersson ( 1992 ) for black holes) to the 
exterior problem with the Lindblom-Detweiler ( Detweiler fc 



Lindblom 1985) scheme for the inside. The key idea is to sep- 



arate ingoing and outgoing solutions by numerically calcu- 
lating their analytic continuations to a place in the complex- 
coordinate plane where they have comparabe amplitudes. 
This is a new approach that could prove to be of great im- 
portance, especially for problems involving rotating stars. In 
such problems the exterior spacetime is not known analyt- 
ically and previous methods (such as the WKB method of 
Kokkotas and S chutz (1992) and the two methods employed 



by Leins et al. (1993)) will consequently fail. But it seems 



likely that a method based on complexifying the coordinates 
will work even if the spacetime itself can only be approxi- 
mated. In fact, the method proposed here may well prove to 
be the only one that remains useful for the strongly damped 
modes of rotating stars. The present approach also has the 
advantage that, although based on numerical integration - 
and in that sense it provides an arbitrarily high numerical 
precision - it solves for quantities that are directly compa- 
rable to those used in the WKB scheme of Kokkotas and 
Schutz. Hence, in this first application of the new method, 
we hope to understand the reason for th e disc repancies be- 
tween the results of Kokkotas and Schutz (1992) and Leins et 
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I. (1993) and thus contribute further to the understanding 



of the oscillation modes of nonrotating relativistic stars. 



2 NONRADIAL PULSATIONS OF A 
RELATIVISTIC STAR 

In Regge- Wheeler gauge the pertur bed metric for a rela- 
tivist ic stellar model takes the form (Detweiler & Lindblom 
198^) 



ds 



2iLur'^'^^ Hie^'^^Ytmdtdr + 
(l - r'Hoe'^'Yim) dr^ + 
{l-/Ke'^'Ye^) (dd'' + si: 



(1) 



where Ho, Hi and K are functions of r only, since we have 
assumed a harmonic dependence on time. We also have the 
standard definition of M{r): 

which plays the role of an effective mass inside radius r. 

As is well-known, a stellar model in equilibrium is de- 
scribed by the Tolman-Oppenheimer-Volkov equations [see 
f or ex ample equations (A3)-(A5) of Lindblom and Detweiler 
(1983)]. These are a system of three coupled first order dif- 
ferential equations that determine the mass M, the metric 
quantity v and the pressure p. Their solution requires that 
an equation of state p — p(p) for the energy density is spec- 
ified. In the present investigation we will assume a simple 
polytropic equation of state 

p = Kp^ , (3) 

where k — lOOkm'^ and we are using units in which c — G — 
1 throughout this paper. This choice of equation of state 
may not be very realistic, but it simplifies the calculations 
somewhat. It also enables u s to c ompare ou r resu lts directly 
to those of Balbinski et al. ( |198E| ), K ojima ( |l988[ ), Kokkotas 
and Schutz (|l992[) and Leins et al. (|l993D. 



Inside the star, the perturbed fluid is described by a 
Lagrangian displacement where 

e ' We Yem , 



d 



i — 2f r iujt ^ -tr 



:Ve'^'^Yt, 



(4) 
(5) 

(6) 



On general grounds one would think that perturbations 
of a spherical star have four degrees of freedom. Two of these 
are associated with the fluid and the remaining two corre- 
spond to the gravitational waves. This means that the flve 
functions Ho, Hi, K, W and V should not be ind ependent. 
As was shown by Detweiler and Lindblom (1985) a system 
of equations that is free from singularities can be obtained 
by defirnng 



X = Lu'^{p + p)e 
l(p + p)e-'^'Ho 



-I//2 



V 



1 / (i'-A)/2 

p e^ " 



W + 



(7) 



where a prime denotes a derivative with respect to r. (A 



misprint in equation (6) of Detweiler and Lindblom ( 198 



has been corrected here.) Then it follows, as a consequence 
of Einstein's equations, that 

(2M + nr + Q)Ho = STrr^e'"^'^ X - 
[(n + l)Q-c.Ve-('^+^']//i + 



2 3- 

nr — LO r e 



-Q{2M -r + Q) 



K 



where 



= M + A-KT^p , 



and 



(8) 



(9) 



(10) 



The interior problem now reduces to a system of four 
first order differential equations for Hi, K , W a nd X [equa- 
tions (8)-(ll) in Detweiler & Lindblom (|l985|)]. Only two 



of the four linearly independent solutions to this system are 
well-behaved at the centre of the star (at r = 0). Further- 
more, the perturbed pressure must vanish at the surface 
(r — R), which implies that X{R) — 0. These conditions 
specify a single acceptable solution for each frequency to. 
Physically, this solution describes the response of the star 
when gravitational waves of the given frequency are inci- 
dent upon it. 

In the exterior of the star the fiuid perturbations van- 
ish, and the two metric perturbations Hi and K can b e 
combined in such a way [see for example Fackerell (1971)] 
that one obtains a single second-order differential equation 
known as the Zerilli equation. This equation can be written 



Z = , 



(11) 



with the effective potential Vz given by 
2M\ 



Vz(r) = 2 (l - —J X 
ri^in + l)r^ + Sn^Mr^ -|- 9nM^r + 9M^ 



The tortoise coordinate r* is defined by 
2M\ d 



dr, V 



= 1- 



r J dr 



(12) 



(13) 



Clearly, this equation allows two linearly independent solu- 
tions. Far away from the star one of these can be identified 
as an outgoing wave, whereas the other describes an ingo- 
ing wave. For a general frequency, the physically acceptable 
solution for the interior of the star leads to a mixture of 
out- and ingoing waves at spatial infinity. The quasinormal 
modes of the system are those special frequencies for which 
no waves come in from infinity. In a sense, these are the 
frequencies at which the star can be expected to radiate 
"spontaneously" . 



3 THE PHASE-AMPLITUDE APPROACH 

In any numerical approach it is imperative that the quan- 
tities under consideration are smooth and nicely behaving. 
If that is not the case, an integration scheme may require 
ridiculously small steps in order to achieve the desired ac- 
curacy. One standard way to avoid this difficulty when one 
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has to solve for (say) an oscillating function - such as the 
wave-like solutions in the stellar problem - is to appeal to 
general properties of the solutions and express them in terms 
of slowly varying functions. In the following we adopt this 
approach to deal with the exterior problem for stars. 

Although the tortoise coordinate r« appears naturally 
[cf. ( p^l )] in a formal description of the perturbation equa- 
tions it may obscure a numerical solution of the problem, 
especially as we plan to use complex coordinates. In order 
to avoid any difficulties we introduce a new dependent vari- 
able ^ according to 



(14) 



The n w e get a new differential equation (Andersson et al. 
199CI 



+ (7(r) 



* = 0, 



where 



U{r) 



-Vz(r) + 



2M 3M^ 



(15) 



(16) 



One of the important properties of an equation such as 
( p^ is that the Wronskian of any two linearly independent 
solutions must be a constant. This means that the general 
solution can always be constructed as a combination of two 
basic solutions that has the form 



-1/2 

• q exp 



±i / qdr 



(17) 



where the function g(r) is a solution to the nonlinear equa- 
tion 



2q dr'^ 



3 
4g2 



dr 



+q-U^O 



(18) 



At first sight, it may seem as if we have replaced the 
relatively simple second-order differential equation ( p^ ) with 
a more complicated equation. In principle that may be true, 
but it should be quite obvious that the function g - as de- 
termined from (Il8|) - may be slowly varying even when the 
solution to ( p^ oscillates wildly. 

Moreover, it is easy to make a connection to the WKB 
approximation here. Whenever [/ is a slowly varying func- 
tion of r one can neglect the two derivatives in (p^. That 

TJ^^'^ f^inpp t.Viat jg tViP raap fnr Ig 



■an nap q 



we shal l use this approximation to generate initial values for 
the integration of (|l^) far away from the star. This is done 
in such a way that g ^ as r — > -f-cxD, i.e., the function ip^ 
represents an ingoing wave whereas t/j" is outgoing. 



4 STOKES PHENOMENON - WHY THE WKB 
APPROACH BREAKS DOWN 

The purpose of this section is to study to what extent a 
WKB approximation of the Zerilli function at the surface of 
the star can be trusted. Specifically, we want to investigate 
whether a geometrical optics argument that the waves un- 
dergo no reflection in the e xterio r spacetime [as was assumed 
by Kokkotas and Schutz ( 1992 )] can be used in this prob- 
lem. If that is the case, one can simply express the exterior 



solution at the surface in terms of U^^'^{R). In general, one 
would expect this kind of argument to be valid for high fre- 
quencies. In the case of really low frequencies, however, the 
situation becomes more obscure. Then waves can possibly 
be backscattered by the curvature of spacetime. 

These qualitative expectations agree with the results of 
Leins et al. (1993) who found a few highly damped modes 
with a very small real part that Kokkotas and Schutz were 
unable to find. On the other hand, Leins et al. did not con- 
firm the " kink" that occurred for relatively high frequencies 
in the WKB spectra. They suggested that this indicates that 
the WKB approach becomes less reliable as the oscillation 
frequency increases. This seems unreasonable: If the WKB 
argument is at all valid, the approximation should become 
more accurate as the oscillation frequency increases. 

It is clear that the exact phase-amplitude approach out- 
lined in the previous section provides us with the means 
to test the WKB assumption quantitatively. By comparing 
the numerically determined function q to the approximation 
U^^^ we can assess the reliability of the method used by 
Kokkotas and Schutz. (It should be noted that this approx- 
imation is not exactly that used by Kokkotas and Schutz 
( [1992| ) where the Zerilli equation was approached directly. 
Nevertheless, a study such as the present one will provide 
results relevant also for that approximation.) 

In approaching the problem for the stellar exterior nu- 
merically we must first devise a scheme that avoids a prob- 
lem that arises for rapidly damped modes. The desired 
outgoing-wave solution to (|l^ will grow exponentially as 
r -|-oo. This means that it is difficult to resolve the in- 
going solution (that is exponentially small) given a finite 
numerical precision. This problem is well-known from stud- 
ies of black-hole normal modes, and can be avoided by let- 
ting r assume complex values [see Andersson (1992)]. Since 
tp~ ^ exp(— zojr) as r — > -l-oo it is clear that the exponential 
divergence can be suppressed along a path in the complex r- 
plane. For large ||r|| the preferred path is a straight line with 
slope given by — Im uj/Ke ui. Such paths are parallell to the 
so-called anti-Stokes lines that play a cruc ial role in a com- 
plex WKB analysis ( Andersson et al. 1993 ). The anti-Stokes 
lines are curves along which J is a real quantity. In 

order to deal with the present problem we will integrate ( |l8| ) 
along a straight line (with the prescribed slope) from a point 
in the asymptotic regime towards the stellar surface. One 
can show, either using an analytic continuation argument or 



a WK B analysis such as that of Araiijo et al. (Araujo et al 



199J), that asymptotic conditions introduced for complex 
values of r in the way described above are, in fact, identical 
to the desired outgoing-wave boundary condition on the real 
r-axis. 

In a complex-coordinate WKB approach to ( |l5| ) the 
so-called Stokes lines also play an i mportant role [see the 
discussion by Andersson et al. (1993)]. These correspond to 
J U^^^dr being purely imaginary. From each of the, possibly 
complex, zeros of U^^'^ - the transition points of the problem 
- emanate three such contours. When an approximate solu- 
tion to ( p^ is continued across a Stokes line the character of 
the solution changes. This is called the Stokes phenomenon. 
Specifically it means that, if a certain linear combination of 
t/)^ and tp" (with q replaced by U^^'^ ) represents the desired 
solution to (hsl) on one side of the Stokes line, another linear 
combination should probably be used on the opposite side. 
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and Schutz cannot yet be dismissed. However, we can pro- 
ceed one step further and use the numerically determined 
function q in the necessary matching at the stellar surface. 
This should provide an accurate way of dealing with the ex- 
terior problem, and it will hopefully enable us to conclude 
whether the "kink" is real or not. 



Figure 1. Top: Comparing the numerically determined function q 
(solid line) to the WKB approximation U^^^ (dashed line). Only 
the real parts of the manifestly complex functions are shown, 
but the result is similar for the imaginary parts. This example 
is for model 4 (see the beginning of section 6.f for particulars) 
and ujM = 0.35 -I- 0.84i. The WKB approach of Kokkotas and 
Schutz cannot be trusted in this case. The numerical calculation 
towards the surface of the star is performed along a straight line 
at an angle of roughly -67 degrees from the real coordinate axis. 
This corresponds to the asymptotic direction of the so-called anti- 
Stokes lines. The Stokes phenomenon gives rise to oscillations in 
q as one gets close to the surface. Bottom: For comparison we 
show results of a similar study on the real r-axis. It should be 
emphasized that the case shown here is one of the most difficult, 
and that the agreement between the exact and the approximate 
function is generally good when Re uj > Im uj. 



Moreover, it is stra ightforward to show [see the Appendix 
in Andersson (199S)] that this effect is important also for 
numerical integration. If the integration of (1^ is continued 
across a Stokes line the function q may become oscillatory. 

We thus have a simple diagnostic of t he "ri o-reflection" 
assumption used by Kokkotas and Schutz ( |l992|) : If a Stokes 
line (emanating from one of the complex zeros of U^^^) 
crosses the real r-axis between r = R and r — +oo we should 
have reflection. For the stellar models considered here we 
find that this does indeed happen for the new modes discov- 
ered by Leins et al. (1993). For small oscillation frequencies 
a Stokes line crosses the real r-axis in the relevant inter- 
val and the Stokes phenomenon should be accounted for in 
an approximate analysis. It is, of course, possible that this 
would give rise to a very small correction. To test whether 
this is the case, we have compared q as numerically deter- 
mined from ( |l8[ ) to U^^^ for several frequencies. A typical 
example - for one o f the new modes that were identified 
by Leins et al. (1993) - is displayed in Figure 1. From this 
Figure it follows that the WKB method used by Kokkotas 
and Schutz ( 199^ ) cannot be trusted for modes with a very 
small real part and a relatively large imaginary part. 

In order to illustrate the difficulties involved in numeri- 
cal integration along the real r-axis - and thus the advantage 
of using complex coordinates - we also show results of a real- 
r calculation in Figure 1. This specific example corresponds 
to the worst possible situation. The numerically computed 
function q has a singularity close to the real r-axis. That 
this can be expected if one fails to account for the Stokes 
phenomenon is clear from the equations given in the Ap- 
pendix of Andersson (1993). Note also that the magnitude 
of q drops dramatically immediately after passing close to 
this singularity. When the surface of the star is reached \q\ 
has decreased more than ten orders of magnitude. 

Our study shows that the WKB technique should be 
reliable for modes with large real parts. This conclusion 
would mean that the "kink" in the spectrum of Kokkotas 



5 AN ACCURATE CONDITION FOR 
NORMAL MODES 

In principle, it is a simple task to derive a normal-mode 
condition based on the phase-amplitude approach. We want 
to match the numerical solution to the inside problem to 
the solution for the exterior at the surface of the star. That 
is, we require that the Zerilli function and its derivative be 
continuous at the surface R. Then we can use (^) and 

-) * 

which follows immediately from (p^. 

Let us first assume that the physically acceptable solu- 
tion for the interior of the star corresponds to a mixture of 
out- and ingoing waves at spatial infinity. Then we have 



dZ _ / 2Afy/^ d* M / 2M 
dr, \ r J dr \ r 



(19) 



outside the star, and we get from ( [l7| ) 

1 dq 
2qd^ 



dr 



1 dq 
"'^YqTr 



(20) 



(21) 



The next step is to solve (|l^) for a particular value of ui 
from a point far away from the star in the way described in 
the previous section. This means that the above expressions 
for ^ and its derivative in terms of will remain valid. 
Moreover, if we take the lower limit of integration in 

(0) to 

be the surface R of the star (where we do the matching) we 
get 



^±(i?)=g-i/^(i?) 



(22) 



which means that we need only determine q and its deriva- 
tive. 

It should be pointed out that the amplitudes Ai-n and 
Aout as defined by ( [20[ ) differ slightly from those that follow 
from the asymptotic behaviour of the Zerilli function; 



Bin 



+ -Boute 



(23) 



The difference is essentially a phase-factor and is of no im- 
portance for a search for mode-frequencies. With our defi- 
nition the amplitudes of the functions that asymptotically 
represent out- and ingoing waves are equal at the surface of 
the star [see (H)]. 

Straightforward algebra leads to 

-1/2 

Ain = ^= (1-^1 X 



2My 



l_dq 



^1 + TT-r 



M_ 

i?2 



+ Z's 



(24) 



and 



A, 



2^/q{R) 



1 - 



2M\ 



-1/2 
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Zs 



/ 2M\ / I dq\ 



+ Z's 



(25) 



where Zs represents the interior solution at the surface of 
the star, and a prime denotes a derivative with respect to 
r*. However, since the solution is only determined up to a 
constant factor by the physical conditions it makes more 
sense to study 



_dq\ _ 
2q drj 



M ' 



Z's 



Zs[{l 



2M ^ 
R J 



(ia + J-^) + M. 

V'i ^ 2q dr) + fl2 



]+Z's 



(26) 



A quasinormal mode of the stellar system corresponds to a 
singularity of this ratio as a function of u. 

It may be worthwhile to clarify the meaning of dq/dr in 
the above equations, since we are assuming that r is com- 
plex in the exterior. In principle, one could simply interpret 
(letting z represent the complex r) dq/dr as dq/dz but this 
may not be very practical. As mentioned previously, we in- 
tegrate for the exterior solution along a straight line in the 
complex z-plane. This line is parametrized by the real dis- 
tance p (from the surface of the star) with constant phase 
angle 6. Then dz/dp = e*^ and dq/dp is the derivative of q 
with respect to the real parameter. Thus we get 

dq ^ dq ^ ^-ie dq 
dr dz dp 

which is easily implemented numerically. 



(27) 



Table 1 . A sample of characteristic frequencies for stellar modes 
associated with model 4 of Kokkotas and Schutz (1992). The 
highly damped (w) modes are compared to the slowly damped 
(p) modes. The results for w-modes are in good agreement with 
those of Leins et al. (1993) The frequencies found for the f-mode 
and the first of the p-modes are in perfect agreement with those 
obtained by Kojima (1988). The "new" w-modes found by Leins 
et al. are indicated by f (the imaginary part of the one in paren- 
thesis is too large for our program to provide truly reliable re- 
sults, but there does seem to be a mode there). The same data 
is shown in Figure 2. All entries are given in units of M~^. The 
given damping rates for the slowly damped modes are accurate 
to a few parts in 10^. 
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6 DISCUSSION OF NUMERICAL RESULTS 



The main purpose of the present study is to understand the 
discre pancies between the resul ts of Kokkotas and Schutz 
(1992) and those of Leins et al. (1993). Specifically, we want 
to find the "new" modes that seem to exis t for very small os - 
cillation frequencies and large damping ( Leins et al. 1993 ). 
We want to see if our ad hoc explanation for why those 
modes could not be found when a WKB approach to the ex- 
terior problem was used (see section 4) is correct, i.e., that 
one must account for Stokes phenomenon for frequencies 
close to the imaginary ti;-axis. We also feel that it is impera- 
tive that we check the possib le existence of the "kink" that 
Kokkotas and Schutz (1992) found in their spectra, espe - 
cially since this feature was not found by Leins et al. (1993). 



6.1 Highly damped (w) modes 

We have performed detailed calcula tions for the four models 
studied by Kokkotas and Schutz (1992). A sample of the 



numerical results for model 4 (that has characteristics: pc = 
10^*^ g/cm^, R = 6.465 km and M = l.SMg, i.e., 2M/R = 
0.594) can be found in Table 1. 

We set out on this investigation to verify the existence of 
modes with large damping and small oscillation frequencies 
that had been found by Leins et al. (1993), and also to find 
whether there is a "kink" in the w-mode spectrum or not. 
The first of these questions was readily answered as soon 
as we had combined the numerical integration approach for 
the exterior vacuum (see sections 3-5) with the numerical 
code that Kokkotas and Schutz had use d for the interior 
problem. The modes found by Leins et al. (11993) do, indeed, 



exist. In view of our understanding of the failure of the WKB 
approach for these frequencies this makes sense. As is clear 
from Figure 1, the WKB approach does break down when 
the real part of t he fr equency is small. 

Leins et al. ( 1995 ) assume d tha t the "kink" in the spec- 
tra of Kokkotas and Schutz (1992) was evidence that the 
WKB method fails also as Re lj increases. We have found 
that this is not the case. Our first calculations gave results 
in good agreement with those of Kokkotas and Schutz for 
high frequencies, with a pronounced "kink' in all spectra, 
and it was clear that the WKB method is reliable for high 
oscillation frequencies. Leins et al. ( 1993| ) argue that their 
approach to the exterior problem did not involve approx- 
imations and should therefore be considered as more reli- 
able than the WKB approximation. Still, we found that our 
new numerical scheme gave results that did not agree well 
with those of Leins et al. At the same time, both meth- 
ods have proved to be reliable and providing results of high 
accuracy in the case of black holes. Hence, it seemed likely 
that they should both be trusted, and that any discrepan- 
cies, eg., the "kink", were due to differences in the approach 
to the interior problem. Unfortunately, this meant that we 
faced a rather d ifficu lt situation. The approach of Detweiler 
and Lindblom ( 1985 ) has been used in all studies of highly 



damped oscillations of relativistic stars. Hence, any differ- 
ences would be in the numerical codes and rather difficult 
to find. 

With this in mind we scrutinized our code for the in- 
terior problem, and found points where it could be im- 
proved. The amended code is much more reliable, and nu- 
merically accurate, than the original one used by Kokkotas 
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and Schutz. It also runs considerably faster. It turns out 
that most of the improvements have a minor effect on the 
numerical results, however. (The results quoted in Table 1 
were, obviously, computed with the final version.) After test- 
ing the code we have a much clearer idea of its restrictions. 
Most importantly, we have found that the "kink" is inti- 
mately related to the choice of point close to the centre 
of the star where power series expansions are used to ini- 
tiate the numerical integration. When choosing a starting 
point one introduces an artificial length-scale in the prob- 
lem, and this scale is reflected by the oscillation frequency 
at which the "kink" occurs. This can easily be illustrated 
if one considers the eigenfunctions: As the oscillation fre- 
quency increases the number of nodes in each eigenfunction 
increases. One will eventually reach a frequency where nodes 
should occur for smaller radii than the one where the numer- 
ical integration is initiated. The failure of the power-series 
solutions to account for this feature of the physical solu- 
tions leads to the "kink" in the spectrum. By initiating the 
integration sufficiently close to the centre of the star one 
can ensure that this "kink" does not occur in the frequency 
regime of interest. In the calculations discussed here we ini- 
tiated the integration at a point roughly co rresp onding to 
2.5% of R. It seems likely that Leins et al. (1993) initiated 



their integration at a point closer to the centre than we did 
originally, and had they continued their calculation to higher 
oscillation frequencies they would h ave fo und a "kink" sim- 
ilar to that of Kokkotas and Schutz (1992). The "kinks" are 
numerical artefacts. 

We also find that it is extremely difficult to avoid nu- 
merical noise in the eigenfunctions. This is mainly due to 
the fact that the eigenfunction X decreases rapidly towards 
the surface of the star (where it should vanish). When solv- 
ing for the interior of the star one matches a solution from 
the centre (that satisfies the physical condition of regular- 
ity) to one from the surface of the star at some intermediate 
point. In order to get a reliable result it is imperative that 
this matching is not done at too large a radius. We find 
that the interior matching should be performed for radii 
smaller than something like 0.25i?. Furthermore, the calcu- 
lation gets more difficult as Im u) increases. This is due to 
the fact that ||X|| is then many magnitudes smaller than 
any of the other three eigenfunctions K, H\ and W . Thus it 
is difflcult to compute X with acceptable accuracy. For this 
reason the calculation is more difficult for less relativistic 
stellar models: The damping of the w- modes generally in- 
creases as the model gets less relativistic. Interesti ngly, this 
explains why the results of Kokkotas and Schutz (1992) for 



model 1 were not as good as for the other three models. In 
view of this fact we decided to restrict our detailed study 
to models 2-4 of Kokkotas and Schutz (1992), and also not 
to attempt a search for modes with Im ujM > 1. We feel 
that the interior problem must be reformulated if a reliable 
search for modes with very large damping is to be made. In 
this context a study based on t he alte rnative formulations 
of Ch andrasekhar and Ferrari (1991a), or Ipser and Price 



(1991), may provide interesting results. However, this does 



not mean that our results should not be taken seriously. We 
believe that the present study is the most accurate one to 
date within the range of validity that we have indicated. 



Figure 2. The p- and w-mode spectra for model 4, cf. Table 1. 
Note the absence of a "kink" for high oscillation frequencies in 
the w-mode spectrum, and that the damping does not decrease 
monotonically as one proceeds up the p-mode spectrum. 



6.2 Slowly damped (p) modes 

An interesting effect of the improvements that we made to 
our code for the interior problem is that it is now consid- 
erably more accurate for modes with very slow damping. 
Hence, it is not at all difficult to compute the f-mode for 
the various models. Furthermore, we are able to iterate for 
several of the first p-modes. As far as we know noone has 
attempted to compute these modes before, and the results 
that we obtain are somewhat surprising. The damping of 
the p-modes does not necessarily, as has been generally pre- 
sumed, decrease rapidly as one proceeds to higher oscillation 
frequencies. In fact, for the most compact of our models 
(model 4) we find an intriguing behaviour. There is a low- 
order mode that damps away much faster than expected (see 
Table 1). This more rapidly damped mode occurs only for 
the most compact of our stellar models. 

This is an interesting result that could be of potential 
importance for many of our astrophysical expectations, such 
as the life-time of gravitational-wave sources. Furthermore, 
methods used by other authors to study f-modes could be 
used to test its correctness. For model 4 the ratio of real 
to imaginary part for the relevant mode is not considerably 
different from that for the f-mode. Therefore, it seems likel y 
that the program used b y Detweiler and Lindblom (1985), 
or that of Kojima (1988), could be used to test our results. 
In fact, Dr Kojima and Dr Lindblom have both very kindly 
performed this calculation. The results they obtain for the 
p-modes of model 4 are in good agreement with those listed 
in our Table 1. 

It may be worthwhile to comment on the difficulty of 
finding p-modes here. It is clear that a numerical code must 
be very robust and accurate if iteration for these slowly 
damped modes is to be at all possible. The imaginary part 
is typically at least six orders of magnitude smaller than 
the real part, and if one requires single-digit precision in the 
imaginary part one must therefore achieve at least six digits 
in the real part. Our amended code should be reliable for 
imaginary parts larger than lO"''. 

In cases when the imaginary part of the mode-frequency 
is so small that iteration is not possible, one can always 
attempt to infer the real part from a graphical approach. 
We know that the modes correspond to the zeros of the 
asymptotic amplitude Ai-n (or the singularities of the ra- 
tio AovLt/Ai-n). If such a zero is situated close to the real 
oj-axis it should, in principle, be easy to distinguish in a 
plot of log||Ain||. This idea wa s rece ntly exploited by, for 
example, Ferrari and Germano ( 199^ ). That it works nicely 
is illustrated in Figure 3. When generating this figure we 
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Figure 3. Illustration of the graphic approach for finding p- 
modes. We show logjg ||j4in|| as a function of real ujM for model 
4. A slowly damped mode corresponds to a narrow singularity, 
and all the p-modes listed in Table 1 can easily be distinguished. 
Peaks in this figure indicate the existence of the w-modes. The 
real parts of the w and p-modes are shown as vertical lines at the 
top and bottom of the figure, respectively. 



normalised the perturbed quantities in such a way that the 
radial displacement W is equal to 1 at the surface of the 
star for all frequencies. The slowly damped p-modes give 
rise to very narrow singularities that are easy to distinguish 
although the data in the figure correspond to AuM — 10""^. 
It is not very time consuming to create a figure with this 
kind of resolution, and it is certainly worthwhile. In fact, 
we can not only pin down all the p-modes given in Table 
1, but also guess at which oscillation frequencies there will 
be a w-mode. It seems as if a peak in log ||Ain|| is an indi- 
cation of a w-mode higher up in the complex t<;-plane. We 
did not expect to find this feature in this kind of graph and 
the exact origin of these peaks is something of a mystery. 
Interestingly the peaks do not show up in similar graphs 
generated by the numerical code of Lindblom (private com- 
munication). It thus seems that the difference between the 
two approaches, i.e. between (i) matching the interior solu- 
tion to the exterior one at the surface of the star and (ii) 
using the value of the Zerilli function at the surface as ini- 
tial data for integration towards infinity where the solution 
is matched to an asymptotic form, is crucial. We believe 
that the information contained in this sort of graph is im- 
portant to our understanding of the stellar pulsations and 
are presently investigating this issue further. Although we do 
not yet have the a clearcut explanation of the phenomenon 
it is worthwhile to make one final point. It is easy to show 
that the peaks in Figure 3 originate from the interior solu- 
tion and not the exterior one. Remembering that the Zerilli 
function is generated from the two spacetime perturbations 
K and Hi, this agrees well with the notion that the space- 
time perturbations dominate the fluid ones (which set the 
normalisation for the figure) at the w-mode frequencies. 



6.3 Are there different families of highly damped 
modes? 



Leins et al. ( 1993 ) argued that they had found a new family 
of highly damped modes. Given the present evidence we have 
to assess whether that is the case or not. We have shown why 
Kokkotas and Schutz (1992) could not find modes with high 



damping and small oscillation frequencies with the WKB 
scheme they used. But surely, that in itself does not imply 
that any modes they could not find belongs to a new family? 
Different families of modes must be distinguished for clear 
physical reasons. The difference between the two families (p- 
and g-modes) of slowly damped modes is well established. 



and the w-modes are certainly distinct. But given the knowl- 
edge available at the moment it is not at all clear that a sim- 
ilar split of the highly damped modes into different classes 
makes sense. 

Leins et al. ( |l993| ) put forward two arguments for why 
their new modes belong to a different family. The first relates 
to the number of nodes of the eigenfunctions corresponding 
to the various modes. They find that the number of nodes 
increase systematically as one progresses up the w-mode se- 
quence, but that the situation is not that clear for the new 
modes. In doing this they claim to be studying the "am- 
plitude" of the eigenfunctions. What they actually study is 
the number of nodes in the real part of each eigenfunction. 
An argument solely based on the real part of a manifestly 
complex eigenfunction clearly does not make much sense. 
From, for example, (|l|) it is evident that both the real and 
the imaginary part will be relevant. Hence, the nodes in the 
imaginary part of the eigenfunction must also be studied. 
When that is done the conclusion regarding the higher oscil- 
lation frequencies remains the same, but for the new modes 
the situation is somewhat clearer. The imaginary part of 
each eigenfunction may have a node even if the real part 
does not, and vice versa. However, although it may provide 
a useful diagnostic, the number of nodes in the eigenfunc- 
tions is not a very good measure for distinguishing different 
families of mo des. It is well known [see for example chap- 
ter 17 in Cox (198C)] that different p-modes may have an 
identical number of nodes in the eigenfunctions, and there 
is no reason why that could not be the case also for highly 
damped modes. 



However, Leins et at (1993) argue that their new mode s 
are different from the ones of Kokkotas and Schutz (1995) 
because of other features of the eigenfunctions. For the new 
modes the eigenfunctions die out towards the centre of the 
star, whereas the amplitude stay roughly constant for all 
values of r for the old modes. This difference is similar to 



that between g- and p-modes (Cox 198C). For the g-modes 
- which correspond to oscillation frequencies smaller than 
that of the f-mode - the eigenfunctions are known to die 
out towards the centre of the star, whereas this is not the 
case for the p-mode eigenfunctions. Could it be that the 
perturbations generally die out towards the centre of the 
star for very small oscillation frequencies? If so, a distinction 
of different w-mode families based on this feature would not 
be satisfactory. 

At the present time, it seems premature to divide the w- 
modes of relativistic stars into different "families" . If such a 
classification is to make sense it must be based on clear phys- 
ical principles, and our understanding of the highly damped 
stellar oscillations is still far from satisfactory. The issue 
could probably be resolved by fur ther detailed stu dies of 
different families of stellar models (Leins et al. 1993). A vi- 
tal piece of information that is still missing from the puzzle 
regards the possible existence of modes with very large imag- 
inary parts. As mentioned previously we have found that 
the present approach to the interior problem suffers from 
numerical difficulties when Im uiM > 1 or so. This is not 
to be taken as an indication that modes with larger imagi- 
nary parts do not exist. On the contrary, such modes may 
well exist and a mode survey that covers the entire complex 
frequency plane could produce interesting results. 
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7 CONCLUDING REMARKS 

In order to understand the discrep ancies between results ob- 
tained by Kokko tas and Schutz (1992) and Leins, Nollert 
and Soffel (1993), we have computed highly damped modes 



for relativistic stars. The stellar models considered are sim- 
ple polytropes used in several previous investigations. We 
have outlined a new numerical approach to the problem for 
the exterior of the star. This new scheme is based on numeri- 
cal integration for complex coordinates, and is similar to one 
t hat has proved extremely reliable for black-hole problems 
( lAndersson 1992|) . 

With this new approach to the stellar exterior we are 
able to find m odes that were not identified by Kokkotas and 
Schutz (1992). These modes, that are highly damped and sit- 
uated close to the imaginary freque ncy- ax is, agree perfectly 
with modes found by Leins et al. ( |l993| ). That the WKB 
method employed by Kokkotas and Schutz failed to distin- 
guish these modes can be understood if the so-called Stokes 
phenomenon (familiar from WKB theory) is accounted for. 

We have managed to explain the occurence of a "kink" 
in the spectra of Kokkotas and Schutz: It is a numerical 
artefact due to the choice of point close to the centre of the 
star where power-series expansions are used to initiate the 
numerical integration. This conclusion was drawn after we 
had scrutinized the way that we dealt with the stellar in- 
terior. We also found other points where our original code 
could be improved. Although none of those changes affected 
the numerical results significantly, the program is now more 
robust and runs considerably faster. That it allows us to do 
accurate calculations is illustrated by the fact that we could 
iterate for several of the very slowly damped p-modes. This 
led to a rather surprising result: For the most compact of 
our stellar models (corresponding to 2M/R = 0.594) the 
damping of the p-modes does not decrease monotonically 
as one proceeds to higher oscillation frequencies. There is a 
low-order mode that damps away at least ten times faster 
than anticipated. This result may significantly affect our ex- 
pectations regarding, for example, stellar stability and the 
lifetime of gravitational-wave sources. Hence, it is of some 
importance that it be studied in more detail. 

The new complex-coordinate method that we have em- 
ployed in the present study seems promising also for rotat- 
ing stars, i.e., when the exterior spacetime is only known 
approximately. Since all previously suggested methods will 
surely fail to handle that difficulty this new approach could 
turn out to be of considerable importance. The anomalous 
p-mode damping that we have discovered here would be es- 
pecially interesting if it exists for rotating star f-modes, be- 
cause these limit stability. Such issues are, however, beyond 
the scope of the present investigation and we will return to 
them in the future. 



provided us with independent checks of our results for slowly 
damped modes. 
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